#Library Setup
library(readr)
library(here)
library(tidyverse)Modeling the Impact of a Household Carbon Tax in Argentina: A Microsimulation Exercise
Loading Dataset
#### Household Income
hh_2018 <- read_delim(here("data","engho2018_hogares.txt"),
delim = "|",
col_names = FALSE,
locale = locale(encoding = "Latin1")
)
colnames(hh_2018) <- strsplit(
"id.provincia.region.subregion.trimestre.anio.pondera.cv1c04.cv1c05_a.cv1c05_b.cv1c05_c.cv1c05_d.cv1c06_e.cv1c06_f.cv1c07.cv1c08.cv1c09.cv1c10.cv1c11.cv1c12.cv1c13.ch01.ch02.ch03.ch04.ch05.ch06.ch07.ch08.ch09.ch10.ch11.ch12.ch13.ch14.ch15.ch16.ch17_a.ch17_b.ch18.ch19.ch20.ch21.propauto.regten.hacina.jniveled.jsexo.jedad_agrup.jsitconyugal.jdif_lp.jestado.jocupengh.jpercept.jcomed.tipohog.cantping.reldep.cantmiem.cantadequi.menor18.menor14.mayor65.dif_lph.cant_act.permieoc.clima_educativo.gastot.gastotpc.gasto_comunh.gascomp.gasvent.gc_01.gc_02.gc_03.gc_04.gc_05.gc_06.gc_07.gc_08.gc_09.gc_10.gc_11.gc_12.gc09_01.gc09_02.gc09_03.gc09_04.gc09_05.gc09_06.gc09_07.gc09_08.gc09_09.fp_contado.fp_credito.fp_especie.fp_prodprop.fp_tarjetas.fp_otras.fp_indef.tn_hipsup.tn_autoserv.tn_negesp.tn_restaurant.tn_kiosco.tn_internet.tn_comedor.tn_especie.tn_otros.tn_indef.decgaphp.decgaphr.decgapht.quigaphp.quigaphr.quigapht.gasto_imputado.m_gasto_imputado.ingtoth.ingpch.m_ingtoth.ingtoth_imp.qinth_p.qinth_r.qinth_t.dinth_p.dinth_r.dinth_t.qinpch_p.qinpch_r.qinpch_t.dinpch_p.dinpch_r.dinpch_t",
"\\."
)[[1]]
# Count of missing values in each column
#colSums(is.na(hh_2018))
# First rows as colnames
colnames(hh_2018) <- as.character(unlist(hh_2018[1, ]))
hh_2018 <- hh_2018[-1, ]
hh_2018 <- type.convert(hh_2018, as.is = TRUE)
##### Household Expenditure
gastos <- read_delim((here("data",
"engho2018_gastos.txt")),
delim = "|",
locale = locale(encoding = "Latin1"))
##### Codebook
codes <- read_delim(here("data",
"engho2018_articulos.txt"),
delim = "|",
locale = locale(encoding = "Latin1")
)Cleaning & Wrangling Data
########################################################
##### Energy expenditure: electricity, gas, others #####
########################################################
energy_exp <- gastos %>%
filter(articulo %in% c("A0451101", "A0451102", "A0451103",
"A0452101", "A0452102", "A0452103",
"A0452104", "A0452105",
"A0453101", "A0453102", "A0453103")) %>%
mutate(articulo_label = case_when(
articulo == "A0451101" ~ "Electricidad - vivienda principal",
articulo == "A0451102" ~ "Electricidad - vivienda secundaria",
articulo == "A0451103" ~ "Gastos de electricidad realizados para una vivienda ocupada por otro hogar",
articulo == "A0452101" ~ "Gas envasado en garrafas",
articulo == "A0452102" ~ "Gas a granel",
articulo == "A0452103" ~ "Gas natural - vivienda principal",
articulo == "A0452104" ~ "Gas envasado en tubo",
articulo == "A0452105" ~ "Gas natural - vivienda secundaria",
articulo == "A0453101" ~ "Kerosene",
articulo == "A0453102" ~ "Leña y carbón",
articulo == "A0453103" ~ "Otros combustibles para el hogar",
TRUE ~ "Otro"
))
###########################
##### Household Income ####
###########################
income_hh <- hh_2018 %>% select(id, ingtoth, ingpch)
## Merging household income with expenditure in energy
carbon_tax <- energy_exp %>%
left_join(income_hh %>% select(id, ingtoth), by = "id") %>%
select (-miembro, -subclase, -clase, -grupo, -division, -forma_pago, -modo_adq, -lugar_adq, -tipo_negocio)subregion_labels <- c(
"1A" = "CABA",
"1B" = "GBA",
"2A" = "Córdoba y La Pampa",
"2B" = "Santa Fe y Entre Ríos",
"2C" = "Buenos Aires (resto)",
"3A" = "Jujuy, Salta y Tucumán",
"3B" = "La Rioja, Catamarca y Santiago del Estero",
"4A" = "Misiones y Corrientes",
"4B" = "Chaco y Formosa",
"5" = "Cuyo (San Juan, Mendoza y San Luis)",
"6A" = "Neuquén y Río Negro",
"6B" = "Chubut, Santa Cruz y Tierra del Fuego"
)Introduction
Motivation
With 394.76 million tonnes of CO₂ emitted in 2021, Argentina accounts for approximately 0.83% of global emissions (Climate Watch 2025). Although not a major polluter, the government implemented a tax reform in 2017 that included the introduction of a carbon tax through the amendment of Law No. 23.966. As a result, Argentina is one of the few Latin American countries with a carbon tax, alongside Colombia, Chile, Mexico, and Uruguay. This reform aligns with a growing global trend toward carbon taxation.
Carbon taxes are increasingly seen as a crucial tool to drive environmental sustainability by reducing greenhouse gas emissions. There are indeed, a cost-effective means of reducing emissions compared to command-and-control regulations or subsidies (Gugler, Haxhimusa, and Liebensteiner 2021). As a Pigouvian tax, carbon pricing helps internalize the external costs of carbon emissions, creating incentives for households to adjust their energy consumption behavior toward cleaner alternatives.
However, literature suggests that they can be regressive, disproportionately burdening lower-income households. This is particularly concerning given the global push toward ecological transition and sustainability. As countries implement carbon tax policies to meet their climate goals, it is essential to assess their potential socioeconomic impacts. Understanding how various income groups will be affected is critical to ensuring that the transition to a sustainable economy is both fair and inclusive. A well-designed policy can help mitigate the risk of exacerbating inequality, ensuring that environmental goals are achieved without disproportionately harming vulnerable populations.
Objective
To evaluate the distributional impact and policy design trade-offs of implementing a household-level carbon tax in Argentina, using a microsimulation model based on household energy expenditure data and carbon pricing structures from Argentina and added policies from developed countries.
Context
Carbon taxes are difficult to design, as they involve the complex task of putting a price on emissions and valuing different sources of carbon dioxide emissions differently. There are many ways to approach this, which is why carbon taxes are regulated in a variety of ways across countries. For the purposes of our microsimulation, we will consider three different carbon tax models, drawing primarily from two countries. We will take inspiration from the Argentinian model, as well as from Canada, which has implemented a more aggressive carbon tax system.
In 2017, Macri’s government introduced a carbon tax, layered on top of an existing tax on liquid fuels. The carbon tax consists of fixed lump sums per unit volume of each fuel, adjusted for inflation every quarter and varying according to the carbon dioxide (CO₂) emissions associated with each fossil fuel. The government planned a gradual increase in the tax rate, an annual 10% hike from 2019 to 2028 (MÉXICO2 - Plataforma Mexicana de Carbono 2023). As of 2023, the carbon price resulting from these increases corresponded to approximately USD 10 per ton of CO₂ (USD 10/tCO₂) (MÉXICO2 - Plataforma Mexicana de Carbono 2023). Tax rates are relatively low compared to the what is advised by the IMF, which sets the floor for middle-income countries at USD 50/tCO₂ (Latin America and Caribbean 2025), but higher than in regional peers such as Chile, Colombia, and Mexico. The tax is applied upstream, targeting fuel importers and producers, rather than directly to households. Its design resembles that of a typical Pigouvian tax.
On the other hand, Canada has a dual carbon pricing system, consisting of a federal fuel charge and a separate system for large industrial emitters. The federal fuel charge is applied to producers, importers, and distributors of fossil fuels, as well as other CO₂-emitting entities such as air and marine carriers (Canada Revenue Agency 2022)[.](https://www.canada.ca/en/revenue-agency/services/tax/excise-taxes-duties-levies/fuel-charge/registration.html)Although the charge is levied upstream, it is often passed on to consumers.
As of April 2025, the federal fuel charge has been discontinued, but its structure will inform the design of this microsimulation (Environment and Climate Change 2018). Like Argentina, Canada sets specific charges for specific fuel types (Canada Department of Finance 2025). However, each province is allowed to implement a higher carbon price than the federal benchmark. Other countries, such as Mexico, also apply differentiated carbon taxes across regions.
Compared to Argentina, Canada enforces a higher and more uniform carbon price (approximately $50/tCO₂e). Additionally, with the exception of the Yukon and Nunavut territories, Canada provides direct lump-sum rebates to households through the Canada Carbon Rebate scheme. These rebates are adjusted based on family size, rural residence, and province, helping to offset the increase in living costs resulting from the carbon tax. This model offers a valuable benchmark for balancing equity and efficiency in carbon pricing policy design.
Moreover, energy sources are not uniform across regions in Argentina. By 2017–2018—the latest available data and the year Argentina introduced a minor carbon tax—an energy subsidy system was in place based on household income and geographic location (province). Subsidies applied to both electricity and natural gas.
Under the 2018 segmentation scheme, residential users were divided into three categories:
- N1 (High income): Pay the full cost of energy.
- N2 (Low income): Receive full subsidies (tarifa social).
- N3 (Middle income): Receive partial subsidies with a consumption cap.
How Final Energy Prices Are Determined Energy bills in Argentina reflect four main components:
Energy generation price – set nationally and uniform across the country. Transmission costs – for transporting energy to distributors. Distribution costs (VAD) – set by local distributors and vary by province/municipality. Taxes – national, provincial, and municipal. Thus, the final price a household pays depends on:
Income level and subsidy eligibility, and The specific distribution and tax rates set by their province or municipality. Electricity vs. Gas Regulation Electricity distribution is decentralized and regulated by provinces and municipalities, except in the Buenos Aires Metropolitan Area (AMBA), where it’s under national jurisdiction. Gas distribution remains under national regulation, managed by ENARGAS, with prices differentiated by service zones and subzones. While the generation cost is the same nationwide, the final bill varies significantly due to differences in local distribution charges (VAD) and regional tax burdens.
Microsimulation Model Structure
Basic Microsimulation Model
Our basic model will have the following characteristics:
Use microdata from Argentina’s national household expenditure survey (Encuesta Nacional de Gastos de los Hogares, ENGHo) to estimate baseline household spending on energy sources (INDEC: Instituto Nacional de Estadística y Censos de la República Argentina 2018).
Apply carbon price increases to each energy source using the official carbon content conversion table (from Argentina’s and Canada’s tax authority).
Model price pass-through to consumers under the assumption that households ultimately bear the tax burden through higher energy prices.mi
As such we have created three basic model that try to emulate three different carbon taxation systems:
A first model that simulates a uniform tax across all households and no inclusion of behavioral response or rebates.
A second model that attempts to mimic the Argentinean carbon tax as it is regulated with lump sums instead of tax rates, also with no behavioral response or rebates.
A third model that includes variation in tax rates depending on how pollutant the energy source used is mimicking the standards from Canada. In this third model we include variations in tax rates by energy type, behavioral response through elasticity, and the inclusion of different schemes of rebates.
Data Description
In the data frame we have constructed basing ourselves on the 2018 household survey, we have information covering 21,065 households, identified by the variable “id”. The data was collected between November of 2017 and November of 2018.
Our data includes spending per household, as well as quantity consumed per household, on the following energy sources:
Bottled gas / LPG in canisters (Gas envasado en garrafas)
Bulk-delivered LPG (Gas a granel)
Gas in larger cylinders or tanks (Gas envasado en tubo)
Markatable Natural Gas (Gas Natural)
Electricity (Electricidad)
Kerosene
Other liquid fossil fuels (Otros combustibles)
Coal and fuelwood (Leña y Carbón)
Sensitivity Analysis
Building on the third model, we will perform a sensitivity analysis that tests multiple scenarios using alternative parameter values.
Elasticity Analysis
Based on our third model we will conduct a study on how increases in carbon prices affect household demand for different energy sources. Drawing on the estimates reported by (Zabaloy and Viego 2022) and choose different levels of elasticities ranging from -0.09 to -0.76. Their findings highlight Argentina’s fuel demand as highly price-inelastic—largely a consequence of long-standing price controls and subsidies, limited short-run substitutes, and behavioural lock-in, with vehicles and appliances calibrated to persistently cheap energy.
Different Carbon Price Levels
We will apply different tax rates to different energy types.
Distributional Analysis
Furthermore, we will evaluate the impact of carbon taxes as set forth on our three initial model on income distribution. In order to do that will first estimate the post-tax energy expenditures for each household, and then estimate the impact on the tax income and income distribution per income quintiles.
We will evaluate revenue recycling mechanisms through cash transfers:
- Flat per capita rebate (uniform transfer to all)
- Progressive rebate (Enhanced rebates for lower-income households.
- Compare net distributional impacts: tax paid minus rebate received.
Microsimulations
Model I: A uniform tax across all households
Model
#Suppose we add a 10% tax on the price paid by all households. What would be the final amount paid, the income effect and revenue effect?
# Build function to calculate carbon tax of 10%
calculate_carbon_tax <- function(monto) {
tax <- monto * 0.10
return(tax)
}
# Apply function to dataframe
carbon_tax_I <- carbon_tax %>%
mutate(
carbon_tax = map_dbl(
.x = monto, #apply the function to every value of the column "monto"
.f = calculate_carbon_tax
),
# Calcular el ingreso neto después del impuesto
net_income = ingtoth - carbon_tax
)
# Calculate quintiles
carbon_tax_I <- carbon_tax_I %>%
# Creamos una variable de decil basada en el ingreso total
mutate(
quintile = ntile(ingtoth, 5)
)
# Avg values by quintile
quintile_averages <- carbon_tax_I %>%
group_by(quintile) %>%
summarise(
avg_total_income = mean(ingtoth, na.rm = TRUE),
avg_net_income = mean(net_income, na.rm = TRUE),
avg_tax = mean(carbon_tax, na.rm = TRUE),
# Calcular el porcentaje que representa el impuesto del ingreso total
tax_percentage = mean(carbon_tax/ingtoth * 100, na.rm = TRUE)
)
# Pivot long for plotting
plot_data <- quintile_averages %>%
select(quintile, avg_total_income, avg_net_income) %>%
pivot_longer(
cols = c(avg_total_income, avg_net_income),
names_to = "income_type",
values_to = "value"
)Output
### OUTPUTS
# $Income – $energy spend: How much income increases or decreases before and after by quintile.
income_plot_model1 <- ggplot(plot_data, aes(x = factor(quintile), y = value, fill = income_type)) +
geom_bar(stat = "identity", position = position_dodge()) +
scale_fill_manual(
values = c("avg_total_income" = "steelblue", "avg_net_income" = "tomato2"),
labels = c("avg_total_income" = "Total Income", "avg_net_income" = "Net Income")
) +
labs(
title = "Effect of Carbon Tax on Household Income by Quintile",
subtitle = "Comparison of average total vs. post-tax (net) income across income groups",
x = "Income Quintile (from lowest to highest)",
y = "Average Monthly Income (in local currency)",
fill = ""
) +
theme_minimal()
print(income_plot_model1)# Relative Tax Burden
tax_burden_plot_model1 <- ggplot(quintile_averages, aes(x = factor(quintile), y = tax_percentage)) +
geom_bar(stat = "identity", fill = "darkgreen") +
labs(
title = "Carbon Tax as a Share of Income by Quintile",
subtitle = "Lower-income households bear a larger burden relative to their income",
x = "Income Quintile (from lowest to highest)",
y = "Carbon Tax (% of Income)"
) +
theme_minimal()
print(tax_burden_plot_model1)#Revenue collected by the government
tax_by_quintile <- carbon_tax_I %>%
group_by(quintile) %>%
summarise(
total_tax = sum(carbon_tax, na.rm = TRUE)
) %>%
ungroup() %>%
mutate(
total_revenue = sum(total_tax),
tax_share = total_tax / total_revenue * 100 # proporción pagada por cada quintil
)
#Contribution by quintile
tax_share_plot_model1 <- ggplot(tax_by_quintile, aes(x = factor(quintile), y = tax_share)) +
geom_bar(stat = "identity", fill = "#5A3989") +
labs(
title = "Carbon Tax Burden by Income Quintile",
subtitle = "Share of government revenue from carbon tax contributed by each income group",
x = "Income Quintile (from lowest to highest)",
y = "Share of Total Carbon Tax Revenue (%)"
) +
theme_minimal()
print(tax_share_plot_model1)## Geographic Analysis
subregion_quintile <- carbon_tax_I %>%
group_by(subregion, quintile) %>%
summarise(
avg_tax_pct = mean(carbon_tax / ingtoth * 100, na.rm = TRUE),
.groups = "drop"
)
subregion_model1 <- ggplot(subregion_quintile, aes(x = factor(quintile), y = subregion, fill = avg_tax_pct)) +
geom_tile() +
scale_y_discrete(labels = subregion_labels) +
scale_fill_viridis_c() +
labs(title = "Carbon Tax Effect by Subregion and Income Quintile",
x = "Income Quintile",
y = "",
fill = "Carbon Tax as a Share of Income") +
theme_minimal()
print(subregion_model1)Model II: Tax varies according to energy source, with more pollutant sources with higher taxes, and is consistent with lump sums established by the Argentinian Government.
In this second version of the model, we aim to create a microsimulation that more closely reflects the carbon taxation system currently in place in Argentina.
Since our household survey is from 2018, with data collected in 2017, we will use the lump sums from the first trimester of 2018 as the basis for our calculations. However, because the survey spans from November 2017 to November 2018—a period during which tax rates changed—we will assume that the applicable lump sum is the average of the four quarterly amounts charged throughout 2018.
Under Argentina’s carbon tax, only a limited number of fossil fuels are taxed. These include:
Unleaded naphtha, up to 92 RON (per litre)
Unleaded naphtha, more than 92 RON (per litre)
Virgin naphtha (per litre) Petroleum solvent (per litre)
White spirit (per litre) Gas oil (per litre)
Diesel oil (per litre) Kerosene (per litre)
Petroleum coke (per kilogram)
Coal (per kilogram)
Among the energy sources used at the household level, only two fall under this taxation scheme: kerosene and the broad category “other liquid fossil fuels.” For these two, we have data on the quantities (in litres) consumed by each household.
Accordingly, we will apply the average of the four quarterly taxes for 2018 to kerosene. For “other liquid fossil fuels,” I will construct an average tax based on the mean of the lump sums applied to all other relevant fuel types. All other energy sources will be considered untaxed under this simulation.
This approach is admittedly imperfect, as it does not account for the share of electricity generated from fossil fuels, notably, coal, which should, in principle, also be subject to carbon taxation. As we don’t know the amounts of these energy sources used to generated that a given amount of electricity, we cannot include them in this model.
As a result, a portion of households will remain untaxed in the model, even though their actual energy consumption would indirectly fall under the carbon tax in real life.
Model
# Build function to calculate differentiated carbon tax rates based on the energy source
calculate_carbon_tax <- function(cantidad, articulo) {
lump_sums <- case_when(
#Only Kerosene and Ohter liquid fossil fuels
articulo == "A0453101" ~ mean(0.473, 0.505, 0.549, 0.626), # Kerosene
articulo == "A0453103" ~ mean( 0.412, 0.440, 0.478, 0.545, 0.412, 0.440, 0.478, 0.545,0.412, 0.440, 0.478, 0.545,0.412, 0.440, 0.478, 0.545,0.412, 0.440, 0.478, 0.545,0.412, 0.505, 0.549, 0.545, 0.473, 0.505, 0.549, 0.626, 0.473, 0.505, 0.549, 0.626), #Otros Combustibles
TRUE ~ 0
)
# Calculate tax
tax <- cantidad * lump_sums
return(tax)
}# Step 1: calculate tax by energy type
carbon_tax_II <- carbon_tax %>%
mutate(
carbon_tax_item = map2_dbl(
.x = cantidad,
.y = articulo,
.f = calculate_carbon_tax
),
# original price paid
price = monto/cantidad,
# per unit markup after tax
tax_per_unit = carbon_tax_item/cantidad,
# Effective tax rate
tax_rate = tax_per_unit / price
)# Step 2: Total tax paid by household
household_tax_II <- carbon_tax_II %>%
group_by(id) %>%
summarise(
total_carbon_tax = sum(carbon_tax_item, na.rm = TRUE),
ingtoth = first(ingtoth) # Take the first income by id (it is the same for each id)
) %>%
mutate(
net_income = ingtoth - total_carbon_tax
)# Step 3: Calculate quintiles based on income
household_tax_II <- household_tax_II %>%
mutate(
quintile = ntile(ingtoth, 5)
)# Step 4: Analyzing income
quintile_analysis_II <- household_tax_II %>%
group_by(quintile) %>%
summarise(
avg_total_income = mean(ingtoth, na.rm = TRUE),
avg_net_income = mean(net_income, na.rm = TRUE),
avg_tax = mean(total_carbon_tax, na.rm = TRUE),
avg_tax_burden = mean(total_carbon_tax/ingtoth * 100, na.rm = TRUE)
)
carbon_tax_with_quintile_II <- carbon_tax_II %>%
left_join(
household_tax_II %>% select(id, quintile, total_carbon_tax),
by = "id"
)Output
#### OUTPUTS
# $Income – $energy spend: How much income increases or decreases before and after by quintile.
income_comparison <- quintile_analysis_II %>%
select(quintile, avg_total_income, avg_net_income) %>%
pivot_longer(
cols = c(avg_total_income, avg_net_income),
names_to = "income_type",
values_to = "amount"
)
income_plot_model2 <- ggplot(income_comparison, aes(x = factor(quintile), y = amount, fill = income_type)) +
geom_bar(stat = "identity", position = position_dodge()) +
scale_fill_manual(
values = c("avg_total_income" = "steelblue", "avg_net_income" = "tomato2"),
labels = c("avg_total_income" = "Total Income", "avg_net_income" = "Net Income")
) +
labs(
title = "Effect of Carbon Tax on Household Income by Quintile",
subtitle = "Comparison of average total vs. post-tax (net) income across income groups",
x = "Income Quintile (from lowest to highest)",
y = "Average Monthly Income (in local currency)",
fill = ""
) +
theme_minimal()
print(income_plot_model2)# Relative Tax Burden
tax_burden_plot_model2 <- ggplot(
quintile_analysis_II,
aes(x = factor(quintile), y = avg_tax_burden)
) +
geom_bar(stat = "identity", fill = "darkgreen") +
labs(
title = "Carbon Tax as a Share of Income by Quintile",
subtitle = "Lower-income households bear a larger burden relative to their income",
x = "Income Quintile (from lowest to highest)",
y = "Carbon Tax (% of Income)"
) +
theme_minimal()
print(tax_burden_plot_model2)# Contribution to gov income revenue
total_revenue_II <- sum(household_tax_II$total_carbon_tax, na.rm = TRUE)
quintile_contributions_II <- household_tax_II %>%
group_by(quintile) %>%
summarise(
total_tax_paid = sum(total_carbon_tax, na.rm = TRUE)
) %>%
mutate(
percent_contribution = 100 * total_tax_paid / total_revenue_II
)
carbon_tax_burden_model2 <- ggplot(quintile_contributions_II, aes(x = factor(quintile), y = percent_contribution)) +
geom_bar(stat = "identity", fill = "#5A3989") +
labs(
title = "Carbon Tax Burden by Income Quintile",
subtitle = "Share of government revenue from carbon tax contributed by each income group",
x = "Income Quintile (from lowest to highest)",
y = "Share of Total Carbon Tax Revenue (%)"
) +
theme_minimal()
print(carbon_tax_burden_model2)#Geographic Analysis
# Step 1: Bring subregion info to household level
household_with_subregion_II <- carbon_tax_II %>%
select(id, subregion) %>%
distinct() %>%
left_join(household_tax_II, by = "id")
# Step 2: Group by subregion and quintile
subregion_quintile_model2 <- household_with_subregion_II %>%
group_by(subregion, quintile) %>%
summarise(
avg_tax_pct = mean(total_carbon_tax / ingtoth * 100, na.rm = TRUE),
.groups = "drop"
)
# Step 3: Plot
subregion_model2 <- ggplot(subregion_quintile_model2, aes(x = factor(quintile), y = subregion, fill = avg_tax_pct)) +
geom_tile() +
scale_y_discrete(labels = subregion_labels) +
scale_fill_viridis_c() +
labs(
title = "Carbon Tax Effect by Subregion and Income Quintile",
x = "Income Quintile",
y = "",
fill = "Carbon Tax as a Share of Income"
) +
theme_minimal()
print(subregion_model2)Model III: Tax varies according to energy source, with more pollutant sources with higher taxes, following the Canadian Model
Description
In this third model, we aim to recreate a more plausible carbon tax system, drawing inspiration from Canada’s approach. Like Argentina, Canada sets absolute rates in Canadian dollars per specific volume or weight of various fossil fuels, with more polluting fuels taxed at higher rates than less polluting ones. It establishes specific rates for the following fossil fuels (Canada Department of Finance 2025):
Aviation gasoline
Aviation turbo fuel
Butane
Coke
Coke oven gas
Combustible waste
Ethane
Gas liquids
Gasoline
Heavy fuel oil
Kerosene
High heat value coal
Kerosene
Light fuel oil
Low heat value coal
Marketable natural gas
Methanol
Naphtha
Non-marketable natural gas
Petroleum coke
Pentanes plus
Propane
Still gas.
As in Argentina, fuel charge rates in Canada have been gradually increasing. The first carbon price introduced in Canada was set at CAD 20/CO₂t in 2019. Rates increased by CAD 10/CO₂t annually until 2022, and by CAD 15/CO₂t each year from 2023 to 2030. Thus, in 2024, the rate reached CAD 75/CO₂t, which is approximately USD 55.55/CO₂t. This is in line with the level at which middle-income countries should be pricing carbon, according to the IMF (USD 50/CO₂t)(Environment and Climate Change 2021).
For the purposes of this study, we make Canada’s fuel charges comparable to those of Argentina. Using the inflation calculator provided by the Bank of Canada (Bank of Canada 2025), we calculate a price increase of 20.24% between 2018 and 2024. Accordingly, we adjust the 2024 fuel charge rates to their 2018 equivalents. We then apply the average 2018 exchange rate between the Argentine peso and the Canadian dollar (1 Argentine Peso = 0.04856 Canadian Dollar (Exchange Rates 2025) to convert the rates into Argentine pesos. It should noted that the exchange rate signficantly fluctuate during this year.
Since not all of the energy sources mentioned above are included in this study, we will group the fossil fuels into broader categories based on the energy sources covered in this survey:
Bottled gas / LPG in canisters (Gas envasado en garrafas): Propane, Butane, Pentanes plus
Bulk-delivered LPG (Gas a granel):Propane, Butane Pentanes plus
Gas in larger cylinders or tanks (Gas envasado en tubo): Propane, Butane, Pentanes plus
Markatable Natural Gas (Gas natural): Natural Gas
Kerosene
Other liquid fossil fuels: Gasoline, Heavy fuel oil, Light fuel oil, Methanol. Naphta
Coal: High heat value coal, Low heat value coal, Coke, Coke oven gas. Since these fuels are rarely used by households bur rather in industrial processes, and the charges associated with them are abnormally high. In this case we are using a tax per kg of 5.6 ARS (0.3 USD at 2018 Argentina’s exchange rate) that reflects a carbon price of roughly USD 125 per tonne CO2 matching the high end of the Paris-aligned benchmark and mirroring countries like Sweden.
As we did in the case of Argentina, we will calculate an average rate per category, corresponding to the mean of the fuel charge rates applied to the relevant energy source usually used by households. For example, the tax applicable to bottled gas will be equal to the average of the rates for Propane, Butane, and Pentanes Plus. Once we have the average rates, we will apply the conversions described above to obtain values that can be compared to Argentina’s 2018 rates.
Code
#Creating fuel charge rates for the different categories of energy sources including in the data by averaging Canada's 2024-2025 fuel charges rates
Canada_rates <- data.frame(
energy_source= c("Propane", "Butane", "Pentanes plus", "Natural Gas", "Kerosene","Gasoline", "Heavy fuel oil", "Light fuel oil", "Methanol", "Naphta"),
rates= c(0.1238, 0.1424, 0.1424, 0.1525, 0.2065, 0.1761, 0.2550, 0.2139, 0.0878, 0.1803))%>%
mutate(
rates = rates/1.2024
)%>%
mutate(
rates = rates/0.04856
)
Canada_rates energy_source rates
1 Propane 2.120279
2 Butane 2.438835
3 Pentanes plus 2.438835
4 Natural Gas 2.611814
5 Kerosene 3.536653
6 Gasoline 3.016003
7 Heavy fuel oil 4.367295
8 Light fuel oil 3.663390
9 Methanol 1.503720
10 Naphta 3.087935
#Averaging these rates into the different energy sources included in the survey
Canada_rate_per_Energy_source_in_survey <- Canada_rates%>%
filter(energy_source == "Propane")%>%
mutate (
Gas_1_Garrafas = mean(c(2.120279, 2.438835, 2.438835 )),
Gas_2_Granel = mean(c(2.120279, 2.438835, 2.438835 )),
Gas_3_Tubo = mean(c(2.120279, 2.438835, 2.438835 )),
Gas_4_Natural = 2.611814,
Kerosene = 3.536653,
Liquid_Fuels = mean(c(3.016003, 4.367295, 3.663390, 1.503720, 3.087935))
)%>%
select(-rates, -energy_source)%>%
pivot_longer(cols = everything(),
names_to = "energy_type",
values_to = "rate")Scenario A
Model A is our pure “status-quo” benchmark: it applies Argentina’s existing carbon levies exactly as they stand—no additional surcharges, no rebates, and no allowance for behavioral response (ε = 0). By holding all tax‐rate increments at zero and assuming perfectly inelastic demand, Model A isolates the mechanical effect of the current tax schedule on household spending shares without any adjustment in consumption or revenue recycling.
Model
elasticity <- 0 # Change for different scenarios
rebate_mode <- "none" # Options: "none", "uniform", "by_quintile"
rebate_target_quintiles <- c(1, 2) # Only used if rebate_mode == "by_quintile"
gas_increment <- 0 # Change ratio to increase / decrease tax on gas articles
fossil_increment <- 0 # Change ratio to increase / decrease tax on fossil articles
firewood_increment <- 0
rebate_share <- 0 # Ratio of revenue that goes into rebate (can't be more than 1)
# Function that applies differentiated tax rates by type of item/energy
calculate_carbon_tax <- function(cantidad, articulo,
fossil = fossil_increment, #ratio of increment
gas = gas_increment,
wood = firewood_increment) #ratio of increment
{
# Tax rates based on type of energy
tax_rates <- case_when(
# Highly polluting fossil fuels
articulo == "A0453101" ~ 3.54* (1 + fossil), # Kerosene
articulo == "A0453103" ~ 3.02* (1 + fossil), # Other fuels
# Gas (medium emissions)
articulo == "A0452101" ~ 2.12 * (1 + gas), # Bottled gas in cylinders
articulo == "A0452102" ~ 2.12 * (1 + gas), # Bulk gas
articulo == "A0452103" ~ 2.61 * (1 + gas), # Natural gas - primary residence
articulo == "A0452104" ~ 2.12 * (1 + gas), # Gas in tubes
articulo == "A0452105" ~ 2.61 * (1 + gas), # Natural gas - secondary residence
# Electricity (varies by energy matrix, but directly less polluting)
articulo == "A0451101" ~ 0, # Electricity - primary residence
articulo == "A0451102" ~ 0, # Electricity - secondary residence
articulo == "A0451103" ~ 0, # Electricity expenses for another residence
# Firewood and coal (may be renewable or not depending on the source)
articulo == "A0453102" ~ 5.6 * (1 + wood), # Firewood and coal
# Others
TRUE ~ 0 # Standard rate
)
# Calculate tax
tax <- cantidad * tax_rates
return(tax)
}# Step 1: calculate tax by energy type
carbon_tax_III_A <- carbon_tax %>%
mutate(
carbon_tax_item = map2_dbl(
.x = cantidad,
.y = articulo,
.f = calculate_carbon_tax
),
#original price paid
price = monto/cantidad,
#per unit markup after tax
tax_per_unit = carbon_tax_item/cantidad,
# Effective tax rate
tax_rate = tax_per_unit / price,
# Apply elasticity adjustment to quantity
cantidad2 = cantidad * (1 + tax_rate) ^ elasticity,
price2 = price * (1 + tax_rate),
monto2 = cantidad2 * price2,
carbon_tax_item2 = cantidad2 * tax_per_unit
)# Step 2: Total tax paid by household
household_tax_III_A <- carbon_tax_III_A %>%
group_by(id) %>%
summarise(
total_carbon_tax = sum(carbon_tax_item, na.rm = TRUE),
total_carbon_tax2 = sum(carbon_tax_item2, na.rm = TRUE),
ingtoth = first(ingtoth), # Take the first income by id (it is the same for each id)
total_spent_energy = sum(monto, na.rm = TRUE),
total_spent_energy2 = sum(monto2, na.rm = TRUE)
) %>%
mutate(
net_income = ingtoth - total_carbon_tax,
net_income2 = ingtoth - total_carbon_tax2
)# Step 3: Calculate quintiles based on income
household_tax_III_A <- household_tax_III_A %>%
mutate(
quintile= ntile(ingtoth, 5)
)
## Step 3.1: Add rebate if any
total_revenue_III_A <- sum(household_tax_III_A$total_carbon_tax, na.rm = TRUE)
total_revenue_III_A2 <- sum(household_tax_III_A$total_carbon_tax2, na.rm = TRUE)
rebate_pool_III_A <- total_revenue_III_A2 * rebate_share
household_tax_III_A <- household_tax_III_A %>%
mutate(
rebate = case_when(
rebate_mode == "uniform" ~ rebate_pool_III_A / n(),
rebate_mode == "by_quintile" & quintile %in% rebate_target_quintiles ~
rebate_pool_III_A / sum(quintile %in% rebate_target_quintiles),
rebate_mode %in% c("by_quintile", "none") ~ 0,
TRUE ~ 0
),
net_income3 = net_income2 + rebate # updated post-rebate income
)# Step 4: Analyzing income
quintile_analysis_III_A <- household_tax_III_A %>%
group_by(quintile) %>%
summarise(
avg_total_income = mean(ingtoth, na.rm = TRUE),
avg_net_income = mean(net_income, na.rm = TRUE),
avg_net_income2 = mean(net_income3, na.rm = TRUE),
avg_tax = mean(total_carbon_tax, na.rm = TRUE),
avg_tax2 = mean(total_carbon_tax2, na.rm = TRUE),
avg_tax_burden = mean(total_carbon_tax/ingtoth * 100, na.rm = TRUE),
avg_tax_burden2 = mean(total_carbon_tax2/ingtoth * 100, na.rm = TRUE)
)
carbon_tax_with_quintile_III_A <- carbon_tax_III_A %>%
left_join(
household_tax_III_A %>% select(id, quintile, total_carbon_tax),
by = "id"
)Outputs
Percent change in energy expenditure as a proportion of income
tax_burden_plot_model_III_A <- household_tax_III_A %>%
group_by(quintile) %>%
summarise(
share_before = sum(total_spent_energy, na.rm = TRUE) / sum(ingtoth, na.rm = TRUE),
share_after = sum(total_spent_energy2, na.rm = TRUE) / sum(ingtoth, na.rm = TRUE)
) %>%
# compute percent change = (after/before – 1) * 100
mutate(
pct_change = (share_after / share_before - 1) * 100
) %>%
ggplot(aes(x = factor(quintile), y = pct_change)) +
geom_bar(stat = "identity", fill = "darkgreen") +
labs(
title = "Percent Change in Energy Spending Share by Quintile",
subtitle = "Post-tax vs. pre-tax share of income (in %)",
x = "Income Quintile (from lowest to highest)",
y = "Change (%)"
) +
theme_minimal()
print(tax_burden_plot_model_III_A)Contribution to gov income revenue
quintile_contributions_III_A <- household_tax_III_A %>%
group_by(quintile) %>%
summarise(
total_tax_paid = sum(total_carbon_tax2, na.rm = TRUE),
total_income = sum(ingtoth),
total_net_income = sum(net_income3)
) %>%
mutate(
percent_contribution = 100 * (total_tax_paid / total_revenue_III_A2),
percent_change_income = 100 * ((total_net_income-total_income)/total_income)
)
carbon_tax_burden_model_III_A <- ggplot(quintile_contributions_III_A, aes(x = factor(quintile), y = percent_contribution)) +
geom_bar(stat = "identity", fill = "#5A3989") +
labs(
title = "Carbon Tax Burden by Income Quintile",
subtitle = "Share of government revenue from carbon tax contributed by each income group",
x = "Income Quintile (from lowest to highest)",
y = "Share of Total Carbon Tax Revenue (%)"
) +
theme_minimal()
print(carbon_tax_burden_model_III_A)Household percent change in net income after tax
net_income_model_III_A <- ggplot(quintile_contributions_III_A, aes(x = factor(quintile), y = percent_change_income)) +
geom_bar(stat = "identity", fill = "#5A3989") +
labs(
title = "Household's net income percent change after tax",
subtitle = "Change in net income by each income group",
x = "Income Quintile (from lowest to highest)",
y = "Percent change in net revenue (%)"
) +
theme_minimal()
print(net_income_model_III_A)Geographic Analysis
# Step 1: Bring subregion info to household level
household_with_subregion_III_A <- carbon_tax_III_A %>%
select(id, subregion) %>%
distinct() %>%
left_join(household_tax_III_A, by = "id")
# Step 2: Group by subregion and quintile
subregion_quintile_model_III_A <- household_with_subregion_III_A %>%
group_by(subregion, quintile) %>%
filter(ingtoth > 0) %>%
summarise(
avg_tax_pct = mean((total_carbon_tax2 / ingtoth) * 100, na.rm = TRUE),
.groups = "drop"
)
# Step 3: Plot
subregion_model_III_A <- ggplot(subregion_quintile_model_III_A, aes(x = factor(quintile), y = subregion, fill = avg_tax_pct)) +
geom_tile() +
scale_y_discrete(labels = subregion_labels) +
scale_fill_viridis_c() +
labs(
title = "Carbon Tax Effect by Subregion and Income Quintile",
x = "Income Quintile",
y = "",
fill = "Carbon Tax as a Share of Income"
) +
theme_minimal()
print(subregion_model_III_A)Scenario B
By imposing a 50% surcharge on fossil fuels, the dirtiest energy sources, kerosene and coal, are hit hardest, while a modest 10% levy on gas only gently nudges users toward cleaner fuel. With a relatively high elasticity (ε = –0.15), this mix generates strong consumption shifts away from taxed fuels and toward subsidized or cleaner alternatives. By forgoing any rebates, the policy keeps the spotlight squarely on pure price signals to drive behavioral change and technology adoption.
Model
elasticity <- -0.15 # Change for different scenarios
rebate_mode <- "none" # Options: "none", "uniform", "by_quintile"
rebate_target_quintiles <- c(1, 2) # Only used if rebate_mode == "by_quintile"
gas_increment <- 0.1 # Change ratio to increase / decrease tax on gas articles
fossil_increment <- 0.5 # Change ratio to increase / decrease tax on fossil articles
firewood_increment <- 0
rebate_share <- 0 # Ratio of revenue that goes into rebate (can't be more than 1)
# Function that applies differentiated tax rates by type of item/energy
calculate_carbon_tax <- function(cantidad, articulo,
fossil = fossil_increment, #ratio of increment
gas = gas_increment,
wood = firewood_increment) #ratio of increment
{
# Tax rates based on type of energy
tax_rates <- case_when(
# Highly polluting fossil fuels
articulo == "A0453101" ~ 3.54* (1 + fossil), # Kerosene
articulo == "A0453103" ~ 3.02* (1 + fossil), # Other fuels
# Gas (medium emissions)
articulo == "A0452101" ~ 2.12 * (1 + gas), # Bottled gas in cylinders
articulo == "A0452102" ~ 2.12 * (1 + gas), # Bulk gas
articulo == "A0452103" ~ 2.61 * (1 + gas), # Natural gas - primary residence
articulo == "A0452104" ~ 2.12 * (1 + gas), # Gas in tubes
articulo == "A0452105" ~ 2.61 * (1 + gas), # Natural gas - secondary residence
# Electricity (varies by energy matrix, but directly less polluting)
articulo == "A0451101" ~ 0, # Electricity - primary residence
articulo == "A0451102" ~ 0, # Electricity - secondary residence
articulo == "A0451103" ~ 0, # Electricity expenses for another residence
# Firewood and coal (may be renewable or not depending on the source)
articulo == "A0453102" ~ 5.6 * (1 + wood), # Firewood and coal
# Others
TRUE ~ 0 # Standard rate
)
# Calculate tax
tax <- cantidad * tax_rates
return(tax)
}# Step 1: calculate tax by energy type
carbon_tax_III_B <- carbon_tax %>%
mutate(
carbon_tax_item = map2_dbl(
.x = cantidad,
.y = articulo,
.f = calculate_carbon_tax
),
#original price paid
price = monto/cantidad,
#per unit markup after tax
tax_per_unit = carbon_tax_item/cantidad,
# Effective tax rate
tax_rate = tax_per_unit / price,
# Apply elasticity adjustment to quantity
cantidad2 = cantidad * (1 + tax_rate) ^ elasticity,
price2 = price * (1 + tax_rate),
monto2 = cantidad2 * price2,
carbon_tax_item2 = cantidad2 * tax_per_unit
)# Step 2: Total tax paid by household
household_tax_III_B <- carbon_tax_III_B %>%
group_by(id) %>%
summarise(
total_carbon_tax = sum(carbon_tax_item, na.rm = TRUE),
total_carbon_tax2 = sum(carbon_tax_item2, na.rm = TRUE),
ingtoth = first(ingtoth), # Take the first income by id (it is the same for each id)
total_spent_energy = sum(monto, na.rm = TRUE),
total_spent_energy2 = sum(monto2, na.rm = TRUE)
) %>%
mutate(
net_income = ingtoth - total_carbon_tax,
net_income2 = ingtoth - total_carbon_tax2
)# Step 3: Calculate quintiles based on income
household_tax_III_B <- household_tax_III_B %>%
mutate(
quintile= ntile(ingtoth, 5)
)
## Step 3.1: Add rebate if any
total_revenue_III_B <- sum(household_tax_III_B$total_carbon_tax, na.rm = TRUE)
total_revenue_III_B2 <- sum(household_tax_III_B$total_carbon_tax2, na.rm = TRUE)
rebate_pool_III_B <- total_revenue_III_B2 * rebate_share
net_revenue_III_B <- total_revenue_III_B2 - rebate_pool_III_B
household_tax_III_B <- household_tax_III_B %>%
mutate(
rebate = case_when(
rebate_mode == "uniform" ~ rebate_pool_III_B / n(),
rebate_mode == "by_quintile" & quintile %in% rebate_target_quintiles ~
rebate_pool_III_B / sum(quintile %in% rebate_target_quintiles),
rebate_mode %in% c("by_quintile", "none") ~ 0,
TRUE ~ 0
),
net_income3 = net_income + rebate # updated post-rebate income
)# Step 4: Analyzing income
quintile_analysis_III_B <- household_tax_III_B %>%
group_by(quintile) %>%
summarise(
avg_total_income = mean(ingtoth, na.rm = TRUE),
avg_net_income = mean(net_income, na.rm = TRUE),
avg_net_income2 = mean(net_income3, na.rm = TRUE),
avg_tax = mean(total_carbon_tax, na.rm = TRUE),
avg_tax2 = mean(total_carbon_tax2, na.rm = TRUE),
avg_tax_burden = mean(total_carbon_tax/ingtoth * 100, na.rm = TRUE),
avg_tax_burden2 = mean(total_carbon_tax2/ingtoth * 100, na.rm = TRUE)
)
carbon_tax_with_quintile_III_B <- carbon_tax_III_B %>%
left_join(
household_tax_III_B %>% select(id, quintile, total_carbon_tax),
by = "id"
)Outputs
Percent change increase in energy expenditure as a proportion of income
pctA <- household_tax_III_A %>%
mutate(
income_after_rebate = ingtoth + rebate
) %>%
group_by(quintile) %>%
summarise(
share_before = sum(total_spent_energy, na.rm = TRUE) / sum(ingtoth, na.rm = TRUE),
share_after = sum(total_spent_energy2, na.rm = TRUE) / sum(income_after_rebate, na.rm = TRUE)
) %>%
mutate(
pct_change = (share_after / share_before - 1) * 100,
pp_change = (share_after - share_before)
)
# 2. Precompute the B pct_change
pctB <- household_tax_III_B %>%
mutate(
income_after_rebate = ingtoth + rebate
) %>%
group_by(quintile) %>%
summarise(
share_before = sum(total_spent_energy, na.rm = TRUE) / sum(ingtoth, na.rm = TRUE),
share_after = sum(total_spent_energy2, na.rm = TRUE) / sum(income_after_rebate, na.rm = TRUE)
) %>%
mutate(
pct_change = ((share_after / share_before) - 1) * 100,
pp_change = (share_after - share_before)
)
pct_change_III_B <- bind_rows(
pctA %>% mutate(model = "A"),
pctB %>% mutate(model = "B")
)
energy_share_change_model_III_B <- ggplot(
pct_change_III_B,
aes(
x = factor(quintile),
y = pct_change,
fill = model
)
) +
geom_bar(
stat = "identity",
position = position_dodge(width = 0.7),
width = 0.6
) +
scale_fill_manual(
values = c(A = "darkgreen", B = "steelblue"),
name = "Model"
) +
labs(
title = "Change in Energy Spending Share by Income Quintile",
subtitle = "Percentage change of energy expenditure as a proportion of income, by model",
x = "Income Quintile (from lowest to highest)",
y = "Change (%)"
) +
scale_x_discrete(labels = as.character(1:5)) +
theme_minimal()
print(energy_share_change_model_III_B)Contribution to gov income revenue
quintile_contributions_III_A <- household_tax_III_A %>%
group_by(quintile) %>%
summarise(
total_tax_paid = sum(total_carbon_tax2, na.rm = TRUE),
total_income = sum(ingtoth, na.rm = TRUE),
total_net_income = sum(net_income3, na.rm = TRUE)
) %>%
mutate(
percent_contribution = 100 * total_tax_paid/ total_revenue_III_A2,
percent_change_income = 100 * ((total_net_income - total_income) / total_income),
model = "A"
) %>%
select(quintile, percent_contribution, percent_change_income, model)
# 2. Summarise B (and tag it)
quintile_contributions_III_B <- household_tax_III_B %>%
group_by(quintile) %>%
summarise(
total_tax_paid = sum(total_carbon_tax2, na.rm = TRUE),
total_income = sum(ingtoth, na.rm = TRUE),
total_net_income = sum(net_income3, na.rm = TRUE),
total_rebate = sum(rebate, na.rm = TRUE)
) %>%
mutate(
percent_contribution = 100 * ((total_tax_paid - total_rebate) / net_revenue_III_B),
percent_change_income = 100 * ((total_net_income - total_income) / total_income),
model = "B"
) %>%
select(quintile, percent_contribution, percent_change_income, model)
quintile_contributions_III_B <- bind_rows(
quintile_contributions_III_A,
quintile_contributions_III_B
)
carbon_tax_burden_model_III_B <- ggplot(
quintile_contributions_III_B,
aes(
x = factor(quintile),
y = percent_contribution,
fill = model
)
) +
geom_bar(
stat = "identity",
position = position_dodge(width = 0.7),
width = 0.6
) +
scale_fill_manual(
values = c(A = "darkgreen", B = "steelblue"),
name = "Model"
) +
labs(
title = "Carbon Tax Burden by Income Quintile",
subtitle = "Share of total carbon tax revenue contributed, by model",
x = "Income Quintile (from lowest to highest)",
y = "Share of Total Carbon Tax Revenue (%)"
) +
scale_x_discrete(labels = as.character(1:5)) +
theme_minimal()
print(carbon_tax_burden_model_III_B)Household percent change in net income after tax
net_income_model_III_B <- ggplot(
quintile_contributions_III_B,
aes(x = factor(quintile),
y = percent_change_income,
fill = model)
) +
geom_bar(
stat = "identity",
position = position_dodge(width = 0.7),
width = 0.6
) +
labs(
title = "Household’s Net Income Change After Tax by Income Quintile",
subtitle = "Percentage change in net income, by model",
x = "Income Quintile (from lowest to highest)",
y = "Percent Change in Net Income (%)",
fill = "Model"
) +
scale_fill_manual(
values = c(A = "darkgreen", B = "steelblue")
) +
theme_minimal()
print(net_income_model_III_B)Geographic Analysis
# Step 1: Bring subregion info to household level
household_with_subregion_III_B <- carbon_tax_III_B %>%
select(id, subregion) %>%
distinct() %>%
left_join(household_tax_III_B, by = "id")
# Step 2: Group by subregion and quintile
subregion_quintile_model_III_B <- household_with_subregion_III_B %>%
group_by(subregion, quintile) %>%
filter(ingtoth > 0) %>%
summarise(
avg_tax_pct = mean((total_carbon_tax2 / ingtoth) * 100, na.rm = TRUE),
.groups = "drop"
)
# Step 3: Plot
subregion_model_III_B <- ggplot(subregion_quintile_model_III_B, aes(x = factor(quintile), y = subregion, fill = avg_tax_pct)) +
geom_tile() +
scale_y_discrete(labels = subregion_labels) +
scale_fill_viridis_c() +
labs(
title = "Carbon Tax Effect by Subregion and Income Quintile",
x = "Income Quintile",
y = "",
fill = "Carbon Tax as a Share of Income"
) +
theme_minimal()
print(subregion_model_III_B)Scenario C
By doubling the fossil‐fuel tax (a 100% surcharge), we send a clear signal to curb the use of the dirtiest fuels—especially kerosene and heavy oils common in small‐scale industry and heating—while a modest 10% gas levy keeps essential cooking fuel affordable for households. With no rebates, we isolate the pure price‐signal effect, maximizing revenue (which can be channeled into clean‐energy investment) and driving significant shifts away from taxed fuels.
Model
elasticity <- -0.05 # Change for different scenarios
rebate_mode <- "none" # Options: "none", "uniform", "by_quintile"
rebate_target_quintiles <- c(1, 2) # Only used if rebate_mode == "by_quintile"
gas_increment <- 0.1 # Change ratio to increase / decrease tax on gas articles
fossil_increment <- 1 # Change ratio to increase / decrease tax on fossil articles
firewood_increment <- 0
rebate_share <- 0 # Ratio of revenue that goes into rebate (can't be more than 1)
# Function that applies differentiated tax rates by type of item/energy
calculate_carbon_tax <- function(cantidad, articulo,
fossil = fossil_increment, #ratio of increment
gas = gas_increment,
wood = firewood_increment) #ratio of increment
{
# Tax rates based on type of energy
tax_rates <- case_when(
# Highly polluting fossil fuels
articulo == "A0453101" ~ 3.54* (1 + fossil), # Kerosene
articulo == "A0453103" ~ 3.02* (1 + fossil), # Other fuels
# Gas (medium emissions)
articulo == "A0452101" ~ 2.12 * (1 + gas), # Bottled gas in cylinders
articulo == "A0452102" ~ 2.12 * (1 + gas), # Bulk gas
articulo == "A0452103" ~ 2.61 * (1 + gas), # Natural gas - primary residence
articulo == "A0452104" ~ 2.12 * (1 + gas), # Gas in tubes
articulo == "A0452105" ~ 2.61 * (1 + gas), # Natural gas - secondary residence
# Electricity (varies by energy matrix, but directly less polluting)
articulo == "A0451101" ~ 0, # Electricity - primary residence
articulo == "A0451102" ~ 0, # Electricity - secondary residence
articulo == "A0451103" ~ 0, # Electricity expenses for another residence
# Firewood and coal (may be renewable or not depending on the source)
articulo == "A0453102" ~ 5.6 * (1 + wood), # Firewood and coal
# Others
TRUE ~ 0 # Standard rate
)
# Calculate tax
tax <- cantidad * tax_rates
return(tax)
}# Step 1: calculate tax by energy type
carbon_tax_III_C <- carbon_tax %>%
mutate(
carbon_tax_item = map2_dbl(
.x = cantidad,
.y = articulo,
.f = calculate_carbon_tax
),
#original price paid
price = monto/cantidad,
#per unit markup after tax
tax_per_unit = carbon_tax_item/cantidad,
# Effective tax rate
tax_rate = tax_per_unit / price,
# Apply elasticity adjustment to quantity
cantidad2 = cantidad * (1 + tax_rate) ^ elasticity,
price2 = price * (1 + tax_rate),
monto2 = cantidad2 * price2,
carbon_tax_item2 = cantidad2 * tax_per_unit
)# Step 2: Total tax paid by household
household_tax_III_C <- carbon_tax_III_C %>%
group_by(id) %>%
summarise(
total_carbon_tax = sum(carbon_tax_item, na.rm = TRUE),
total_carbon_tax2 = sum(carbon_tax_item2, na.rm = TRUE),
ingtoth = first(ingtoth), # Take the first income by id (it is the same for each id)
total_spent_energy = sum(monto, na.rm = TRUE),
total_spent_energy2 = sum(monto2, na.rm = TRUE)
) %>%
mutate(
net_income = ingtoth - total_carbon_tax,
net_income2 = ingtoth - total_carbon_tax2
)# Step 3: Calculate quintiles based on income
household_tax_III_C <- household_tax_III_C %>%
mutate(
quintile= ntile(ingtoth, 5)
)
## Step 3.1: Add rebate if any
total_revenue_III_C <- sum(household_tax_III_C$total_carbon_tax, na.rm = TRUE)
total_revenue_III_C2 <- sum(household_tax_III_C$total_carbon_tax2, na.rm = TRUE)
rebate_pool_III_C <- total_revenue_III_C2 * rebate_share
net_revenue_III_C <- total_revenue_III_C2 - rebate_pool_III_C
household_tax_III_C <- household_tax_III_C %>%
mutate(
rebate = case_when(
rebate_mode == "uniform" ~ rebate_pool_III_C / n(),
rebate_mode == "by_quintile" & quintile %in% rebate_target_quintiles ~
rebate_pool_III_C / sum(quintile %in% rebate_target_quintiles),
rebate_mode %in% c("by_quintile", "none") ~ 0,
TRUE ~ 0
),
net_income3 = net_income + rebate # updated post-rebate income
)# Step 4: Analyzing income
quintile_analysis_III_C <- household_tax_III_C %>%
group_by(quintile) %>%
summarise(
avg_total_income = mean(ingtoth, na.rm = TRUE),
avg_net_income = mean(net_income, na.rm = TRUE),
avg_net_income2 = mean(net_income3, na.rm = TRUE),
avg_tax = mean(total_carbon_tax, na.rm = TRUE),
avg_tax2 = mean(total_carbon_tax2, na.rm = TRUE),
avg_tax_burden = mean(total_carbon_tax/ingtoth * 100, na.rm = TRUE),
avg_tax_burden2 = mean(total_carbon_tax2/ingtoth * 100, na.rm = TRUE)
)
carbon_tax_with_quintile_III_C <- carbon_tax_III_C %>%
left_join(
household_tax_III_C %>% select(id, quintile, total_carbon_tax),
by = "id"
)Outputs
Percent change increase in energy expenditure as a proportion of income
# Precompute the C‐series pct_change
pctC <- household_tax_III_C %>%
mutate(
income_after_rebate = ingtoth + rebate
) %>%
group_by(quintile) %>%
summarise(
share_before = sum(total_spent_energy, na.rm = TRUE) / sum(ingtoth, na.rm = TRUE),
share_after = sum(total_spent_energy2, na.rm = TRUE) / sum(income_after_rebate, na.rm = TRUE)
) %>%
mutate(
pct_change = (share_after / share_before - 1) * 100,
pp_change = (share_after - share_before)
)
pct_change_III_C <- bind_rows(
pct_change_III_B,
pctC %>% mutate(model = "C")
)
energy_share_change_model_III_C <- ggplot(
pct_change_III_C,
aes(
x = factor(quintile),
y = pct_change,
fill = model
)
) +
geom_bar(
stat = "identity",
position = position_dodge(width = 0.7),
width = 0.6
) +
scale_fill_manual(
values = c(A = "darkgreen", B = "steelblue", C = "darkorange"),
name = "Model"
) +
labs(
title = "Change in Energy Spending Share by Income Quintile",
subtitle = "Percentage change of energy expenditure as a proportion of income, by model",
x = "Income Quintile (from lowest to highest)",
y = "Change (%)"
) +
scale_x_discrete(labels = as.character(1:5)) +
theme_minimal()
print(energy_share_change_model_III_C)Contribution to gov income revenue
# Summarise C
quintile_contributions_III_C <- household_tax_III_C %>%
group_by(quintile) %>%
summarise(
total_tax_paid = sum(total_carbon_tax2, na.rm = TRUE),
total_income = sum(ingtoth, na.rm = TRUE),
total_net_income = sum(net_income3, na.rm = TRUE),
total_rebate = sum(rebate, na.rm = TRUE)
) %>%
mutate(
percent_contribution = 100 * ((total_tax_paid - total_rebate) / net_revenue_III_C),
percent_change_income = 100 * ((total_net_income - total_income) / total_income),
model = "C"
) %>%
select(quintile, percent_contribution, percent_change_income, model)
quintile_contributions_III_C <- bind_rows(
quintile_contributions_III_B,
quintile_contributions_III_C
)
carbon_tax_burden_model_III_C <- ggplot(
quintile_contributions_III_C,
aes(
x = factor(quintile),
y = percent_contribution,
fill = model
)
) +
geom_bar(
stat = "identity",
position = position_dodge(width = 0.7),
width = 0.6
) +
scale_fill_manual(
values = c(A = "darkgreen", B = "steelblue", C = "darkorange"),
name = "Model"
) +
labs(
title = "Carbon Tax Burden by Income Quintile",
subtitle = "Share of total carbon tax revenue contributed, by model",
x = "Income Quintile (from lowest to highest)",
y = "Share of Total Carbon Tax Revenue (%)"
) +
scale_x_discrete(labels = as.character(1:5)) +
theme_minimal()
print(carbon_tax_burden_model_III_C)Household percent change in net income after tax
net_income_model_III_C <- ggplot(
quintile_contributions_III_C,
aes(x = factor(quintile),
y = percent_change_income,
fill = model)
) +
geom_bar(
stat = "identity",
position = position_dodge(width = 0.7),
width = 0.6
) +
labs(
title = "Household’s Net Income Change After Tax by Income Quintile",
subtitle = "Percentage change in net income, by model",
x = "Income Quintile (from lowest to highest)",
y = "Percent Change in Net Income (%)",
fill = "Model"
) +
scale_fill_manual(
values = c(A = "darkgreen", B = "steelblue", C = "darkorange")
) +
theme_minimal()
print(net_income_model_III_C)Geographic Analysis
# Step 1: Bring subregion info to household level
household_with_subregion_III_C <- carbon_tax_III_C %>%
select(id, subregion) %>%
distinct() %>%
left_join(household_tax_III_C, by = "id")
# Step 2: Group by subregion and quintile
subregion_quintile_model_III_C <- household_with_subregion_III_C %>%
group_by(subregion, quintile) %>%
filter(ingtoth > 0) %>%
summarise(
avg_tax_pct = mean((total_carbon_tax2 / ingtoth) * 100, na.rm = TRUE),
.groups = "drop"
)
# Step 3: Plot
subregion_model_III_C <- ggplot(subregion_quintile_model_III_C, aes(x = factor(quintile), y = subregion, fill = avg_tax_pct)) +
geom_tile() +
scale_y_discrete(labels = subregion_labels) +
scale_fill_viridis_c() +
labs(
title = "Carbon Tax Effect by Subregion and Income Quintile",
x = "Income Quintile",
y = "",
fill = "Carbon Tax as a Share of Income"
) +
theme_minimal()
print(subregion_model_III_C)Scenario D
In this scenario we introduce modest carbon surcharges, 20% on gas, 15% on fossil fuels and 10% on firewood, while allowing a moderate behavioral response (ε = –0.08). Those price signals will nudge households to trim energy use without triggering drastic cutbacks, and by recycling 25% of all revenues back to consumers via a uniform rebate, we soften the burden equally across the board. The combination of differentiated surcharges and a broad-based rebate tests a revenue-neutral design: it still discourages pollution through higher prices, captures some of the resulting efficiency gains in the form of clean-energy investment funds, and cushions low- and middle-income families from seeing their energy bills spike.
Model
elasticity <- -0.08 # Change for different scenarios
rebate_mode <- "uniform" # Options: "none", "uniform", "by_quintile"
rebate_target_quintiles <- c(1, 2) # Only used if rebate_mode == "by_quintile"
gas_increment <- 0.2 # Change ratio to increase / decrease tax on gas articles
fossil_increment <- 0.15 # Change ratio to increase / decrease tax on fossil articles
firewood_increment <- 0.1
rebate_share <- 0.25 # Ratio of revenue that goes into rebate (can't be more than 1)
# Function that applies differentiated tax rates by type of item/energy
calculate_carbon_tax <- function(cantidad, articulo,
fossil = fossil_increment, #ratio of increment
gas = gas_increment,
wood = firewood_increment) #ratio of increment
{
# Tax rates based on type of energy
tax_rates <- case_when(
# Highly polluting fossil fuels
articulo == "A0453101" ~ 3.54* (1 + fossil), # Kerosene
articulo == "A0453103" ~ 3.02* (1 + fossil), # Other fuels
# Gas (medium emissions)
articulo == "A0452101" ~ 2.12 * (1 + gas), # Bottled gas in cylinders
articulo == "A0452102" ~ 2.12 * (1 + gas), # Bulk gas
articulo == "A0452103" ~ 2.61 * (1 + gas), # Natural gas - primary residence
articulo == "A0452104" ~ 2.12 * (1 + gas), # Gas in tubes
articulo == "A0452105" ~ 2.61 * (1 + gas), # Natural gas - secondary residence
# Electricity (varies by energy matrix, but directly less polluting)
articulo == "A0451101" ~ 0, # Electricity - primary residence
articulo == "A0451102" ~ 0, # Electricity - secondary residence
articulo == "A0451103" ~ 0, # Electricity expenses for another residence
# Firewood and coal (may be renewable or not depending on the source)
articulo == "A0453102" ~ 5.6 * (1 + wood), # Firewood and coal
# Others
TRUE ~ 0 # Standard rate
)
# Calculate tax
tax <- cantidad * tax_rates
return(tax)
}# Step 1: calculate tax by energy type
carbon_tax_III_D <- carbon_tax %>%
mutate(
carbon_tax_item = map2_dbl(
.x = cantidad,
.y = articulo,
.f = calculate_carbon_tax
),
#original price paid
price = monto/cantidad,
#per unit markup after tax
tax_per_unit = carbon_tax_item/cantidad,
# Effective tax rate
tax_rate = tax_per_unit / price,
# Apply elasticity adjustment to quantity
cantidad2 = cantidad * (1 + tax_rate) ^ elasticity,
price2 = price * (1 + tax_rate),
monto2 = cantidad2 * price2,
carbon_tax_item2 = cantidad2 * tax_per_unit
)# Step 2: Total tax paid by household
household_tax_III_D <- carbon_tax_III_D %>%
group_by(id) %>%
summarise(
total_carbon_tax = sum(carbon_tax_item, na.rm = TRUE),
total_carbon_tax2 = sum(carbon_tax_item2, na.rm = TRUE),
ingtoth = first(ingtoth), # Take the first income by id (it is the same for each id)
total_spent_energy = sum(monto, na.rm = TRUE),
total_spent_energy2 = sum(monto2, na.rm = TRUE)
) %>%
mutate(
net_income = ingtoth - total_carbon_tax,
net_income2 = ingtoth - total_carbon_tax2
)# Step 3: Calculate quintiles based on income
household_tax_III_D <- household_tax_III_D %>%
mutate(
quintile= ntile(ingtoth, 5)
)
## Step 3.1: Add rebate if any
total_revenue_III_D <- sum(household_tax_III_D$total_carbon_tax, na.rm = TRUE)
total_revenue_III_D2 <- sum(household_tax_III_D$total_carbon_tax2, na.rm = TRUE)
rebate_pool_III_D <- total_revenue_III_D2 * rebate_share
net_revenue_III_D <- total_revenue_III_D2 - rebate_pool_III_D
household_tax_III_D <- household_tax_III_D %>%
mutate(
rebate = case_when(
rebate_mode == "uniform" ~ rebate_pool_III_D / n(),
rebate_mode == "by_quintile" & quintile %in% rebate_target_quintiles ~
rebate_pool_III_D / sum(quintile %in% rebate_target_quintiles),
rebate_mode %in% c("by_quintile", "none") ~ 0,
TRUE ~ 0
),
net_income3 = net_income + rebate # updated post-rebate income
)# Step 4: Analyzing income
quintile_analysis_III_D <- household_tax_III_D %>%
group_by(quintile) %>%
summarise(
avg_total_income = mean(ingtoth, na.rm = TRUE),
avg_net_income = mean(net_income, na.rm = TRUE),
avg_net_income2 = mean(net_income3, na.rm = TRUE),
avg_tax = mean(total_carbon_tax, na.rm = TRUE),
avg_tax2 = mean(total_carbon_tax2, na.rm = TRUE),
avg_tax_burden = mean(total_carbon_tax/ingtoth * 100, na.rm = TRUE),
avg_tax_burden2 = mean(total_carbon_tax2/ingtoth * 100, na.rm = TRUE)
)
carbon_tax_with_quintile_III_D <- carbon_tax_III_D %>%
left_join(
household_tax_III_D %>% select(id, quintile, total_carbon_tax),
by = "id"
)Outputs
Percent change increase in energy expenditure as a proportion of income
# Precompute the C‐series pct_change
pctD <- household_tax_III_D %>%
mutate(
income_after_rebate = ingtoth + rebate
) %>%
group_by(quintile) %>%
summarise(
share_before = sum(total_spent_energy, na.rm = TRUE) / sum(ingtoth, na.rm = TRUE),
share_after = sum(total_spent_energy2, na.rm = TRUE) / sum(income_after_rebate, na.rm = TRUE)
) %>%
mutate(
pct_change = (share_after / share_before - 1) * 100,
pp_change = (share_after - share_before)
)
pct_change_III_D <- bind_rows(
pct_change_III_C,
pctD %>% mutate(model = "D")
)
energy_share_change_model_III_D <- ggplot(
pct_change_III_D,
aes(
x = factor(quintile),
y = pct_change,
fill = model
)
) +
geom_bar(
stat = "identity",
position = position_dodge(width = 0.7),
width = 0.6
) +
scale_fill_manual(
values = c(A = "darkgreen", B = "steelblue", C = "darkorange", D = "firebrick"),
name = "Model"
) +
labs(
title = "Change in Energy Spending Share by Income Quintile",
subtitle = "Percentage change of energy expenditure as a proportion of income, by model",
x = "Income Quintile (from lowest to highest)",
y = "Change (%)"
) +
scale_x_discrete(labels = as.character(1:5)) +
theme_minimal()
print(energy_share_change_model_III_D)Contribution to gov income revenue
# Summarise D
quintile_contributions_III_D <- household_tax_III_D %>%
group_by(quintile) %>%
summarise(
total_tax_paid = sum(total_carbon_tax2, na.rm = TRUE),
total_income = sum(ingtoth, na.rm = TRUE),
total_net_income = sum(net_income3, na.rm = TRUE),
total_rebate = sum(rebate, na.rm = TRUE)
) %>%
mutate(
percent_contribution = 100 * ((total_tax_paid - total_rebate) / net_revenue_III_D),
percent_change_income = 100 * ((total_net_income - total_income) / total_income),
model = "D"
) %>%
select(quintile, percent_contribution, percent_change_income, model)
quintile_contributions_III_D <- bind_rows(
quintile_contributions_III_C,
quintile_contributions_III_D
)
carbon_tax_burden_model_III_D <- ggplot(
quintile_contributions_III_D,
aes(
x = factor(quintile),
y = percent_contribution,
fill = model
)
) +
geom_bar(
stat = "identity",
position = position_dodge(width = 0.7),
width = 0.6
) +
scale_fill_manual(
values = c(A = "darkgreen", B = "steelblue", C = "darkorange", D = "firebrick"),
name = "Model"
) +
labs(
title = "Carbon Tax Burden by Income Quintile",
subtitle = "Share of total carbon tax revenue contributed, by model",
x = "Income Quintile (from lowest to highest)",
y = "Share of Total Carbon Tax Revenue (%)"
) +
scale_x_discrete(labels = as.character(1:5)) +
theme_minimal()
print(carbon_tax_burden_model_III_D)Household percent change in net income after tax
net_income_model_III_D <- ggplot(
quintile_contributions_III_D,
aes(x = factor(quintile),
y = percent_change_income,
fill = model)
) +
geom_bar(
stat = "identity",
position = position_dodge(width = 0.7),
width = 0.6
) +
labs(
title = "Household’s Net Income Change After Tax by Income Quintile",
subtitle = "Percentage change in net income, by model",
x = "Income Quintile (from lowest to highest)",
y = "Percent Change in Net Income (%)",
fill = "Model"
) +
scale_fill_manual(
values = c(A = "darkgreen", B = "steelblue", C = "darkorange", D = "firebrick")
) +
theme_minimal()
print(net_income_model_III_D)Geographic Analysis
# Step 1: Bring subregion info to household level
household_with_subregion_III_D <- carbon_tax_III_D %>%
select(id, subregion) %>%
distinct() %>%
left_join(household_tax_III_D, by = "id")
# Step 2: Group by subregion and quintile
subregion_quintile_model_III_D <- household_with_subregion_III_D %>%
group_by(subregion, quintile) %>%
filter(ingtoth > 0) %>%
summarise(
avg_tax_pct = mean((total_carbon_tax2 / ingtoth) * 100, na.rm = TRUE),
.groups = "drop"
)
# Step 3: Plot
subregion_model_III_D <- ggplot(subregion_quintile_model_III_D, aes(x = factor(quintile), y = subregion, fill = avg_tax_pct)) +
geom_tile() +
scale_y_discrete(labels = subregion_labels) +
scale_fill_viridis_c() +
labs(
title = "Carbon Tax Effect by Subregion and Income Quintile",
x = "Income Quintile",
y = "",
fill = "Carbon Tax as a Share of Income"
) +
theme_minimal()
print(subregion_model_III_D)Scenario E
A bolder tax increase (+30%) across all polluting fuels to drive stronger carbon cuts, paired with a large targeted rebate (40% of revenues) for the poorest 40% to protect equity. The higher elasticity amplifies both emission reductions and welfare protections.
Model
elasticity <- -0.12 # Change for different scenarios
rebate_mode <- "by_quintile" # Options: "none", "uniform", "by_quintile"
rebate_target_quintiles <- c(1, 2) # Only used if rebate_mode == "by_quintile"
gas_increment <- 0.3 # Change ratio to increase / decrease tax on gas articles
fossil_increment <- 0.3 # Change ratio to increase / decrease tax on fossil articles
firewood_increment <- 0.3
rebate_share <- 0.4 # Ratio of revenue that goes into rebate (can't be more than 1)
# Function that applies differentiated tax rates by type of item/energy
calculate_carbon_tax <- function(cantidad, articulo,
fossil = fossil_increment, #ratio of increment
gas = gas_increment,
wood = firewood_increment) #ratio of increment
{
# Tax rates based on type of energy
tax_rates <- case_when(
# Highly polluting fossil fuels
articulo == "A0453101" ~ 3.54* (1 + fossil), # Kerosene
articulo == "A0453103" ~ 3.02* (1 + fossil), # Other fuels
# Gas (medium emissions)
articulo == "A0452101" ~ 2.12 * (1 + gas), # Bottled gas in cylinders
articulo == "A0452102" ~ 2.12 * (1 + gas), # Bulk gas
articulo == "A0452103" ~ 2.61 * (1 + gas), # Natural gas - primary residence
articulo == "A0452104" ~ 2.12 * (1 + gas), # Gas in tubes
articulo == "A0452105" ~ 2.61 * (1 + gas), # Natural gas - secondary residence
# Electricity (varies by energy matrix, but directly less polluting)
articulo == "A0451101" ~ 0, # Electricity - primary residence
articulo == "A0451102" ~ 0, # Electricity - secondary residence
articulo == "A0451103" ~ 0, # Electricity expenses for another residence
# Firewood and coal (may be renewable or not depending on the source)
articulo == "A0453102" ~ 5.6 * (1 + wood), # Firewood and coal
# Others
TRUE ~ 0 # Standard rate
)
# Calculate tax
tax <- cantidad * tax_rates
return(tax)
}# Step 1: calculate tax by energy type
carbon_tax_III_E <- carbon_tax %>%
mutate(
carbon_tax_item = map2_dbl(
.x = cantidad,
.y = articulo,
.f = calculate_carbon_tax
),
#original price paid
price = monto/cantidad,
#per unit markup after tax
tax_per_unit = carbon_tax_item/cantidad,
# Effective tax rate
tax_rate = tax_per_unit / price,
# Apply elasticity adjustment to quantity
cantidad2 = cantidad * (1 + tax_rate) ^ elasticity,
price2 = price * (1 + tax_rate),
monto2 = cantidad2 * price2,
carbon_tax_item2 = cantidad2 * tax_per_unit
)# Step 2: Total tax paid by household
household_tax_III_E <- carbon_tax_III_E %>%
group_by(id) %>%
summarise(
total_carbon_tax = sum(carbon_tax_item, na.rm = TRUE),
total_carbon_tax2 = sum(carbon_tax_item2, na.rm = TRUE),
ingtoth = first(ingtoth), # Take the first income by id (it is the same for each id)
total_spent_energy = sum(monto, na.rm = TRUE),
total_spent_energy2 = sum(monto2, na.rm = TRUE)
) %>%
mutate(
net_income = ingtoth - total_carbon_tax,
net_income2 = ingtoth - total_carbon_tax2
)# Step 3: Calculate quintiles based on income
household_tax_III_E <- household_tax_III_E %>%
mutate(
quintile= ntile(ingtoth, 5)
)
## Step 3.1: Add rebate if any
total_revenue_III_E <- sum(household_tax_III_E$total_carbon_tax, na.rm = TRUE)
total_revenue_III_E2 <- sum(household_tax_III_E$total_carbon_tax2, na.rm = TRUE)
rebate_pool_III_E <- total_revenue_III_E2 * rebate_share
net_revenue_III_E <- total_revenue_III_E2 - rebate_pool_III_E
household_tax_III_E <- household_tax_III_E %>%
mutate(
rebate = case_when(
rebate_mode == "uniform" ~ rebate_pool_III_E / n(),
rebate_mode == "by_quintile" & quintile %in% rebate_target_quintiles ~
rebate_pool_III_E / sum(quintile %in% rebate_target_quintiles),
rebate_mode %in% c("by_quintile", "none") ~ 0,
TRUE ~ 0
),
net_income3 = net_income + rebate # updated post-rebate income
)# Step 4: Analyzing income
quintile_analysis_III_E <- household_tax_III_E %>%
group_by(quintile) %>%
summarise(
avg_total_income = mean(ingtoth, na.rm = TRUE),
avg_net_income = mean(net_income, na.rm = TRUE),
avg_net_income2 = mean(net_income3, na.rm = TRUE),
avg_tax = mean(total_carbon_tax, na.rm = TRUE),
avg_tax2 = mean(total_carbon_tax2, na.rm = TRUE),
avg_tax_burden = mean(total_carbon_tax/ingtoth * 100, na.rm = TRUE),
avg_tax_burden2 = mean(total_carbon_tax2/ingtoth * 100, na.rm = TRUE)
)
carbon_tax_with_quintile_III_E <- carbon_tax_III_E %>%
left_join(
household_tax_III_E %>% select(id, quintile, total_carbon_tax),
by = "id"
)Outputs
Percent change increase in energy expenditure as a proportion of income
# Precompute the C‐series pct_change
pctE <- household_tax_III_E %>%
mutate(
income_after_rebate = ingtoth + rebate
) %>%
group_by(quintile) %>%
summarise(
share_before = sum(total_spent_energy, na.rm = TRUE) / sum(ingtoth, na.rm = TRUE),
share_after = sum(total_spent_energy2, na.rm = TRUE) / sum(income_after_rebate, na.rm = TRUE)
) %>%
mutate(
pct_change = (share_after / share_before - 1) * 100,
pp_change = (share_after - share_before)
)
pct_change_III_E <- bind_rows(
pct_change_III_D,
pctE %>% mutate(model = "E")
)
energy_share_change_model_III_E <- ggplot(
pct_change_III_E,
aes(
x = factor(quintile),
y = pct_change,
fill = model
)
) +
geom_bar(
stat = "identity",
position = position_dodge(width = 0.7),
width = 0.6
) +
scale_fill_manual(
values = c(A = "darkgreen", B = "steelblue", C = "darkorange", D = "firebrick", E = "darkmagenta"),
name = "Model"
) +
labs(
title = "Change in Energy Spending Share by Income Quintile",
subtitle = "Percentage change of energy expenditure as a proportion of income, by model",
x = "Income Quintile (from lowest to highest)",
y = "Change (%)"
) +
scale_x_discrete(labels = as.character(1:5)) +
theme_minimal()
print(energy_share_change_model_III_E)Contribution to gov income revenue
# Summarise D
quintile_contributions_III_E <- household_tax_III_E %>%
group_by(quintile) %>%
summarise(
total_tax_paid = sum(total_carbon_tax2, na.rm = TRUE),
total_income = sum(ingtoth, na.rm = TRUE),
total_net_income = sum(net_income3, na.rm = TRUE),
total_rebate = sum(rebate, na.rm = TRUE)
) %>%
mutate(
percent_contribution = 100 * ((total_tax_paid - total_rebate) / net_revenue_III_E),
percent_change_income = 100 * ((total_net_income - total_income) / total_income),
model = "E"
) %>%
select(quintile, percent_contribution, percent_change_income, model)
quintile_contributions_III_E <- bind_rows(
quintile_contributions_III_D,
quintile_contributions_III_E
)
carbon_tax_burden_model_III_E <- ggplot(
quintile_contributions_III_E,
aes(
x = factor(quintile),
y = percent_contribution,
fill = model
)
) +
geom_bar(
stat = "identity",
position = position_dodge(width = 0.7),
width = 0.6
) +
scale_fill_manual(
values = c(A = "darkgreen", B = "steelblue", C = "darkorange", D = "firebrick", E = "darkmagenta"),
name = "Model"
) +
labs(
title = "Carbon Tax Burden by Income Quintile",
subtitle = "Share of total carbon tax revenue contributed, by model",
x = "Income Quintile (from lowest to highest)",
y = "Share of Total Carbon Tax Revenue (%)"
) +
scale_x_discrete(labels = as.character(1:5)) +
theme_minimal()
print(carbon_tax_burden_model_III_E)Household percent change in net income after tax
net_income_model_III_E <- ggplot(
quintile_contributions_III_E,
aes(x = factor(quintile),
y = percent_change_income,
fill = model)
) +
geom_bar(
stat = "identity",
position = position_dodge(width = 0.7),
width = 0.6
) +
labs(
title = "Household’s Net Income Change After Tax by Income Quintile",
subtitle = "Percentage change in net income, by model",
x = "Income Quintile (from lowest to highest)",
y = "Percent Change in Net Income (%)",
fill = "Model"
) +
scale_fill_manual(
values = c(A = "darkgreen", B = "steelblue", C = "darkorange", D = "firebrick", E = "darkmagenta")
) +
theme_minimal()
print(net_income_model_III_E)Geographic Analysis
# Step 1: Bring subregion info to household level
household_with_subregion_III_E <- carbon_tax_III_E %>%
select(id, subregion) %>%
distinct() %>%
left_join(household_tax_III_E, by = "id")
# Step 2: Group by subregion and quintile
subregion_quintile_model_III_E <- household_with_subregion_III_E %>%
group_by(subregion, quintile) %>%
filter(ingtoth > 0) %>%
summarise(
avg_tax_pct = mean((total_carbon_tax2 / ingtoth) * 100, na.rm = TRUE),
.groups = "drop"
)
# Step 3: Plot
subregion_model_III_E <- ggplot(subregion_quintile_model_III_E, aes(x = factor(quintile), y = subregion, fill = avg_tax_pct)) +
geom_tile() +
scale_y_discrete(labels = subregion_labels) +
scale_fill_viridis_c() +
labs(
title = "Carbon Tax Effect by Subregion and Income Quintile",
x = "Income Quintile",
y = "",
fill = "Carbon Tax as a Share of Income"
) +
theme_minimal()
print(subregion_model_III_E)